﻿/********************************************************
 *  ██████╗  ██████╗████████╗██╗
 * ██╔════╝ ██╔════╝╚══██╔══╝██║
 * ██║  ███╗██║        ██║   ██║
 * ██║   ██║██║        ██║   ██║
 * ╚██████╔╝╚██████╗   ██║   ███████╗
 *  ╚═════╝  ╚═════╝   ╚═╝   ╚══════╝
 * Geophysical Computational Tools & Library (GCTL)
 *
 * Copyright (c) 2023  Yi Zhang (yizhang-geo@zju.edu.cn)
 *
 * GCTL is distributed under a dual licensing scheme. You can redistribute 
 * it and/or modify it under the terms of the GNU Lesser General Public 
 * License as published by the Free Software Foundation, either version 2 
 * of the License, or (at your option) any later version. You should have 
 * received a copy of the GNU Lesser General Public License along with this 
 * program. If not, see <http://www.gnu.org/licenses/>.
 * 
 * If the terms and conditions of the LGPL v.2. would prevent you from using 
 * the GCTL, please consider the option to obtain a commercial license for a 
 * fee. These licenses are offered by the GCTL's original author. As a rule, 
 * licenses are provided "as-is", unlimited in time for a one time fee. Please 
 * send corresponding requests to: yizhang-geo@zju.edu.cn. Please do not forget 
 * to include some description of your company and the realm of its activities. 
 * Also add information on how to contact you by electronic and paper mail.
 ******************************************************/

#ifndef _GCTL_GLNI_H
#define _GCTL_GLNI_H

#include "../core.h"
#include "../gctl_config.h"

namespace gctl
{
    namespace glni
    {
        class coefficients
        {
        public:
            coefficients();
            coefficients(size_t size);
            virtual ~coefficients();
            
            void initiate(size_t size);
            void clear();

            double node(size_t idx);
            double weight(size_t idx);

        private:
            size_t size_;
            double *nodes_;
            double *weights_;
        };

        /**
         * @brief   Gauss Legendre Numerical Integration in the reference 1D segment
         * 
         * integral range for the reference element: [-1, 1]
         * -1        1
         * |----------|----> u
         * 
         * a class object must be initiated before use. This could be done by passing a desired node 
         * number of the GLNI to the constructor or by calling the initiate() function. To use 
         * a different order the GLNI, one must firstly calling the clear() function before the class 
         * object could be initiated again. To use the class object for integration, the recommended way 
         * is to use the predefined callback interface for customized functions. Once the functions are 
         * defined. It could be passed to the integral() function for evaluation. Also, one could use the 
         * node() and weight() functions to get the Gauss nodes' positions and weights, respectively. 
         * Note that the nodes' position are located wrt. the reference element. Proper scaling must be 
         * applied to map the nodes to a customized range.
         */
        class line
        {
        public:
            line();
            line(size_t size);
            virtual ~line();
            virtual double line_func(double x, double x1, double x2, void *att) = 0;

            void clear();
            size_t size(){return size_;}
            void initiate(size_t size);

            double integral(double low, double hig, void *att = nullptr);

        private:
            coefficients coeff_;
            bool initialized_;
            size_t size_;
        };

        /**
         * @brief   Gauss Legendre Numerical Integration in the 2D quadrilateral
         * 
         * integral range for the reference element: [-1, 1]x[-1, 1]
         * y
         * |
         * [-1,1]    [1,1]
         * |----------|
         * |          |
         * |          |
         * |          |
         * |----------|----> x
         * [-1,-1]   [1,-1]
         * 
         * a class object must be initiated before use. This could be done by passing a desired node 
         * number of the GLNI to the constructor or by calling the initiate() function. To use 
         * a different order the GLNI, one must firstly calling the clear() function before the class 
         * object could be initiated again. To use the class object for integration, the recommended way 
         * is to use the predefined callback interface for customized functions. Once the functions are 
         * defined. It could be passed to the integral() function for evaluation. Since specified transformations 
         * are required to calculate the integral in a 2D quadrilateral area. No direct access to the 
         * nodes and weights of the GLNI are provided here.
         */
        class quadrilateral
        {
        public:
            quadrilateral();
            quadrilateral(size_t size);
            virtual ~quadrilateral();
            virtual double quadrilateral_func(double x, double y, double x1, double x2, double y1, double y2, void *att) = 0;

            void clear();
            int size(){return size_;}
            void initiate(size_t size);

            double integral(double xlow, double xhig, double ylow, double yhig, void *att = nullptr);

        private:
            coefficients coeff_;
            bool initialized_;
            size_t size_;
        };

        /**
         * @brief   Gauss Legendre Numerical Integration in the 2D triangle
         * 
         * integral range for the reference element: [0, 0]x[0, 1]x[0, 1]
         * v
         * |
         * [0, 1]
         * |\
         * |  \
         * |    \
         * |      \
         * |        \
         * |          \
         * |------------\---> u
         * [0, 0]     [1, 0]
         * 
         * a class object must be initiated before use. This could be done by passing a desired node 
         * number of the GLNI to the constructor or by calling the initiate() function. To use 
         * a different order the GLNI, one must firstly calling the clear() function before the class 
         * object could be initiated again. To use the class object for integration, the recommended way 
         * is to use the predefined callback interface for customized functions. Once the functions are 
         * defined. It could be passed to the integral() function for evaluation. Since specified transformations 
         * are required to calculate the integral in a 2D triangular area. No direct access to the 
         * nodes and weights of the GLNI are provided here.
         */
        class triangle
        {
        public:
            triangle();
            triangle(size_t size);
            virtual ~triangle();
            virtual double triangle_func(double x, double y, double x1, double x2, double x3, double y1, double y2, double y3, void *att) = 0;

            void clear();
            int size(){return size_;}
            void initiate(size_t size);

            double integral(double x1, double x2, double x3, double y1, double y2, double y3, void *att = nullptr);

        private:
            coefficients coeff_;
            bool initialized_;
            size_t size_;
        };

        /**
         * @brief   Gauss Legendre Numerical Integration in the 2D triangle
         * 
         * integral range for the reference element: [0, 0]x[0, 1]x[0, 1]
         * v
         * |
         * [0, 1]
         * |\
         * |  \
         * |    \
         * |      \
         * |        \
         * |          \
         * |------------\---> u
         * [0, 0]     [1, 0]
         * 
         * a class object must be initiated before use. This could be done by passing a desired node 
         * number of the GLNI to the constructor or by calling the initiate() function. To use 
         * a different order the GLNI, one must firstly calling the clear() function before the class 
         * object could be initiated again. To use the class object for integration, the recommended way 
         * is to use the predefined callback interface for customized functions. Once the functions are 
         * defined. It could be passed to the integral() function for evaluation. Since specified transformations 
         * are required to calculate the integral in a 2D triangular area. No direct access to the 
         * nodes and weights of the GLNI are provided here.
         */
        class triangle_c
        {
        public:
            triangle_c();
            triangle_c(size_t size);
            virtual ~triangle_c();
            virtual std::complex<double> triangle_c_func(double x, double y, double x1, double x2, double x3, double y1, double y2, double y3, void *att) = 0;

            void clear();
            int size(){return size_;}
            void initiate(size_t size);

            std::complex<double> integral(double x1, double x2, double x3, double y1, double y2, double y3, void *att = nullptr);

        private:
            coefficients coeff_;
            bool initialized_;
            size_t size_;
        };

        /**
         * @brief   Gauss Legendre Numerical Integration in the reference 3D cube
         * 
         * integral range for the reference element: [-1, -1, -1] ~ [1, 1, 1]
         * z
         * |
         * |    ----------- [1,1,1]
         * |  /|         /|
         * | / |        / |
         * |/          /  |
         * |----------|   |
         * |   |      |   |
         * |  /-------|--/
         * | /        | /
         * |----------|--------> x
         * [-1,-1,-1]
         * 
         * a class object must be initiated before use. This could be done by passing a desired node 
         * number of the GLNI to the constructor or by calling the initiate() function. To use 
         * a different order the GLNI, one must firstly calling the clear() function before the class 
         * object could be initiated again. To use the class object for integration, the recommended way 
         * is to use the predefined callback interface for customized functions. Once the functions are 
         * defined. It could be passed to the integral() function for evaluation. Since specified transformations 
         * are required to calculate the integral in a 3D cubic area. No direct access to the 
         * nodes and weights of the GLNI are provided here.
         */
        class cube
        {
        public:
            cube();
            cube(size_t size);
            virtual ~cube();
            virtual double cube_func(double x, double y, double z, double x1, double x2, double y1, double y2, double z1, double z2, void *att) = 0;

            void clear();
            int size(){return size_;}
            void initiate(size_t size);

            double integral(double x1, double x2, double y1, double y2, double z1, double z2, void *att = nullptr);

        private:
            coefficients coeff_;
            bool initialized_;
            size_t size_;
        };

        /**
         * @brief   Gauss Legendre Numerical Integration in the reference tetrahedron
         * 
         * integral range for the reference element: [0, 0, 0]x[1, 0, 0]x[0, 1, 0]x[0, 0, 1]
         * zta
         * |
         * |  [0, 0, 1]
         * |\
         * |  \
         * |    \
         * |      \
         * |        \
         * |          \ [0, 1, 0]
         * |------------\---> eta
         * | [0, 0, 0] /
         * |         /
         * |       /
         * |     /
         * |   /
         * | /
         * |/ [1, 0, 0]
         * |
         * |
         * ksi
         * 
         * a class object must be initiated before use. This could be done by passing a desired node 
         * number of the GLNI to the constructor or by calling the initiate() function. To use 
         * a different order the GLNI, one must firstly calling the clear() function before the class 
         * object could be initiated again. To use the class object for integration, the recommended way 
         * is to use the predefined callback interface for customized functions. Once the functions are 
         * defined. It could be passed to the integral() function for evaluation. Since specified transformations 
         * are required to calculate the integral in a 3D tetrahedral area. No direct access to the 
         * nodes and weights of the GLNI are provided here.
         */
        class tetrahedron
        {
        public:
            tetrahedron();
            tetrahedron(size_t size);
            virtual ~tetrahedron();
            virtual double tetrahedron_func(double ksi, double eta, double zta, 
                double x1, double x2, double x3, double x4, 
                double y1, double y2, double y3, double y4,
                double z1, double z2, double z3, double z4, void *att) = 0;

            void clear();
            int size(){return size_;}
            void initiate(size_t size);

            double integral(double x1, double x2, double x3, double x4, 
                double y1, double y2, double y3, double y4, 
                double z1, double z2, double z3, double z4, 
                void *att = nullptr);

        private:
            coefficients coeff_;
            bool initialized_;
            size_t size_;
        };

        /**
         * @brief	Gauss Legendre Numerical Integration in the reference triangular prism
         * 
         * zta eta
         * |  /      3(0,1,1)
         * | /      /|\
         * |/      / | \
         * |      /  |  \
         * |     /  / \0(0,1,0)
         * |   4|---------| 5
         * |    | /     \ | 
         * |    |/-------\|
         * |   1(0,0,0)   2(1,0,0) 
         * O-----------------> ksi
         * 
         * a class object must be initiated before use. This could be done by passing a desired node 
         * number of the GLNI to the constructor or by calling the initiate() function. To use 
         * a different order the GLNI, one must firstly calling the clear() function before the class 
         * object could be initiated again. To use the class object for integration, the recommended way 
         * is to use the predefined callback interface for customized functions. Once the functions are 
         * defined. It could be passed to the integral() function for evaluation. Since specified transformations 
         * are required to calculate the integral. No direct access to the nodes and weights of the GLNI are provided here.
         */
        class triprism
        {
        public:
            triprism();
            triprism(size_t size);
            virtual ~triprism();
            virtual double triprism_func(double ksi, double eta, double zta, 
                double *xs, double *ys, double *zs, void *att) = 0;

            void clear();
            int size(){return size_;}
            void initiate(size_t size);

            double integral(double *xs, double *ys, double *zs, void *att = nullptr);

        private:
            coefficients coeff_;
            bool initialized_;
            size_t size_;
        };
    }
}

#endif // _GCTL_GLNI_H